Convective instability and mass transport of diffusion layers in a Hele-Shaw geometry 
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We consider experimentally the instability and mass transport of a porous-medium flow in a Hele-Shaw 
geometry. In an initially stable configuration, a lighter fluid (water) is located over a heavier fluid (propylene 
glycol). The fluids mix via diffusion with some regions of the resulting mixture being heavier than either pure 
fluid. Density-driven convection occurs with downward penetrating dense fingers that transport mass much more 
effectively than diffusion alone. We investigate the initial instability and the quasi steady state. The convective 
time and velocity scales, finger width, wave number selection, and normalized mass transport are determined 
for 6, 000 < Ra < 90, 000. The results have important implications for determining the time scales and 
rates of dissolution trapping of carbon dioxide in brine aquifers proposed as possible geologic repositories for 
sequestering carbon dioxide. 



PACS numbers: 47.15.gp, 47.20.Bp, 47.56.+r 

Convective processes that transport heat and mass are criti- 
cal to many natural phenomena in the atmosphere and oceans 
as well as in planetary and stellar interiors H). Another im- 
portant example of current scientific and technological inter- 
est is the convective mixing of supercritical carbon dioxide 
sequestered in a brine-saturated porous medium, [2| i.e., a 
deep saline aquifer. When injected into such aquifers, the 
low density CO2 rises to the top of the aquifer where it is 
contained by cap rock formations. As a pure phase, the buoy- 
ancy of the supercritical CO2 may drive it back to the surface 
defeating attempts to sequester it. Diffusion into the under- 
lying brine results in a denser mixture eliminating buoyancy, 
but this process is very slow and the pure phase of supercriti- 
cal CO2 would remain mostly intact for long times (« 10,000 
years) unless other processes occur. 

One scenario for enhanced dissolution of the CO2 relies on 
this unusual property of the C02-brine system, i.e. C02-brine 
mixtures are denser than pure brineQ making gravitationally- 
driven convection is possible. However, the diffusion layer 
must thicken enough to become unstable, a process similar 
to the thickening and subsequent instability of thermal layers 
||3l . Analytical calculations and numerical simulations have 
been performed |4}j9) to determine the instability incubation 
time and the subsequent mass transport in C02-brine systems. 
Corresponding laboratory experiments lfTOl . however, have not 
succeeded in capturing the essential features of such diffusion- 
driven, convectively unstable layers. 

The difficulties in realizing a laboratory experiment for the 
purposes of testing the scenario of unstable diffusion layers in 
porous media are considerable, e.g. establishing the stable ini- 
tial conditions assumed by calculations and simulations. Also, 
direct experiments using supercritical CO2 and brine are com- 
plicated by high pressure (wlOO bar). Instead, we have devel- 
oped an atmospheric-pressure analog fluid system that mimics 
many features of the CC>2-brine system. It consists of water 
located vertically above propylene glycol (PPG) in a Hele- 



Shaw cell. As an analog system, water replaces the CO2 and 
propylene glycol replaces the brine. In this cell, we can create 
relatively unperturbed initial conditions such that the instabil- 
ity occurs uniformly over a wide lateral portion of the cell. 
We measure the temporal evolution of the system using op- 
tical shadowgraph and determine the wave number selection, 
the time and velocity scales, the finger width, and the mass 
transport rate. Our results establish the basic description of 
diffusively-driven convective instability in such systems. 

The experimental apparatus in Fig. [T] consists of a 
vertically-oriented Hele-Shaw cell with two 1.2-cm-thick 
polycarbonate sheets separated by steel shims of thickness b. 
The flow in the narrow gap b between the parallel surfaces 
obeys Darcy's Law with an isotropic permeability K = b 2 /12 
and a porosity of unity. We vary b from 0.25 to 0.48 mm 
giving values of if in the range 5,000-20,000 Darcy (i.e. 
0.54 x 10~ 4 < K < 1.94 x 10~ 4 cm 2 ). The lateral dimen- 
sion of the flow region is 7.6 cm, and the initial depth of the 
bottom fluid H ranges from 1.25 to 5 cm. Water (p w = 1.000 
g/cc) is lighter than PPG (p vvg = 1.035 g/cc[ 11]) and is ini- 
tially above the PPG. Water-PPG mixtures with a water con- 
centration c w between and w 50% have p > pppa- After a 
period of diffusive growth, the mixed fluid interface may be- 
come gravitationally unstable, similar to the proposed insta- 
bility when supercritical CO2 dissolves into brine. There are, 
however, some important differences between the two fluid 
systems. Water and PPG are fully miscible whereas CO2 is 
partially miscible in brine with a low saturation concentra- 
tion so that the transport properties of the CO2 -brine system 
do not depend strongly on the CO2 concentration. In con- 
trast, the viscosity and mass diffusivity of the water-PPG sys- 
tem are highly dependent on concentration. When comput- 
ing experimental parameters or making comparisons to cal- 
culations or simulations, we use the transport properties of a 
30% water-70% PPG mixture, i.e., v = 0.1 cm 2 /s flU and 
D = 1 x 10~ 6 cm 2 /s (estimated using data from lfT2l ). which 
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FIG. 1: a) Schematic illustration of the Hele-Shaw cell where gap 
width b, PPG reservoir height H, and vertical coordinate z are indi- 
cated. A thin metal shim separating the two fluids is withdrawn hor- 
izontally to begin an experimental run. Images at three intermediate 
times t* are shown for a run with K=1.9 x 10 -4 cm 2 and Ra=87,000: 
b) t* = 2.0, c) t* = 3.0 , d) t* = 10.0. Only aportion of the vertical 
extent of the image is shown. As seen in d), light scattering off of the 
metallic shims makes the precise location of the interface uncertain. 
The vertical arrow in d) indicates the nucleation site of a new finger. 

corresponds to the mixture of maximum density p= 1.044 g/cc 
or the maximum density difference of Ap = 0.009 g/cc relative 
to the underlying PPG. 

The upper and lower sections of the Hele-Shaw cell are sep- 
arated by a system of thin shims of maximum thickness 0.01 
cm. The upper and lower reservoirs are filled with pure wa- 
ter and pure PPG, respectively. A 0.0076-cm-thick "shutter 
shim" is slowly removed horizontally from between the two 
fluids leaving behind the 0.01 -cm-thick "frame shim" to en- 
sure a liquid-tight seal. Small perturbations of the fluid-fluid 
interface were observed during some runs, but these did not 
seem to influence the flow appreciably for measurements pre- 
sented here. The removal of the shutter shim opens up a nar- 
row 0.01 -cm-thick horizontal gap at the interface between the 
upper and lower cells. Our cell filling method ensures that this 
pocket is occupied by water after the shim removal fl3l . 

The development of the flow is visualized using optical 
shadowgraph which images the mixture concentration via its 
index of refraction variation. The images resolve an area of 
lateral size L = 3.0 cm and height 2.4 cm with 1280 x 1024 
pixels in the horizontal and vertical directions, respectively. 
Approximately 2000 frames are gathered during each run with 
the camera frame rate dependent on the permeability of the 
cell. (Total run times range between 2.5 and 10 hours). The 
flow develops through several stages. Initially, the PPG pref- 
erentially diffuses into the water owing to an asymmetric de- 
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FIG. 2: (Color online) a) 5/ X vs z/S. Different values of H are 
H = 1.27 cm (blue), H = 2.54 cm (red), and H = 5.08 cm 
(black). The different symbols within each color indicate differ- 
ent values of K. Over the range of the experiments in this work, 
5 is found to be roughly the same for all runs with an average 
value of 5 — 0.065 cm. The solid line is the best fit of the form 
S/X = (1.90 + 0.038z/<5)/(5.30 + z/<5). b) Area of water consumed 
A/(LH) vs t* for several data sets with K = 0.77 x 10" 4 cm 2 and 
Ra — 8, 700 (open triangles), Ra = 17, 400 (solid squares), and 
Ra = 34, 800 (open circles). The small offset in the initial A/LH 
is due to the uncertainty in the location of the initial fluid-fluid inter- 
face as described in caption of Fig.[T] 

pendence of D(C) 11121 . As the dense region in the diffusion 
layer thickens and the influence of buoyancy overcomes lat- 
eral diffusion, hydrodynamic instability occurs as small per- 
turbations in the dense diffusive interface grow into downward 
moving fingers. Figure [lj>-d show the evolution of a typical 
experimental run at various times after initiation. 

The equations governing the conservation of momentum, 
mass, and concentration; their various non-dimensional scal- 
ings; and solution techniques have been described in detail 
elsewhere [4-9]. Here, we discuss the dimensionless param- 
eters and point out a few details relevant to our experiment. 
Under the conditions of the experiment, the fluids and their 
mixture are incompressible. Conservation of momentum is 
adequately described by a two-dimensional Darcy's law be- 
cause the Reynolds numbers i?e& based on b and the maxi- 
mum observed velocities are less than 10 -3 . The evolution 
of the mixture concentration is governed by the advection- 
diffusion equation, although the molecular diffusion coeffi- 
cient D may be modified by Taylor dispersion [14] in certain 
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regions of higher flow velocity. For example, near the tips and 
within the core of fingers which move downward with velocity 
Vf, the Peclet number Pe — bvf/D m 50 increasing the ef- 
fective lateral diffusion coefficient by a factor 15. Elsewhere, 
the fluid velocities are significantly lower, and the effect on D 
is negligible. 

Recalling that Re^ <C 1 and neglecting the impact of Pe 
within the finger core, the only dimensionless number that 
governs the flow is the Rayleigh number which we define 
in the same way as [4-9], i.e., Ra = ApgKH / pvD where 
g is the acceleration due to gravity and, as described above, 
all fluid properties are evaluated at maximum mixture den- 
sity. Next, we define characteristic length, velocity, and time 
scales. The vertical convective motion of a blob of maximum 
density fluid is given by a balance between buoyancy and vis- 
cous drag, i.e., v c = ApgK/ pv. At v c , the time required for 
this blob to fall the entire depth H of the lower fluid in Fig. 
[TJis t c = Hjv c = H/{ApgK/ pv). We define t* = t/t c . 
Comparative time and velocity scales are the diffusion time 
td = H 2 /D, i.e., the approximate time to establish a diffusion 
gradient over H, and the diffusion velocity Vd = D/H. From 
these scaling relations, we see that Ra = td/t c = v c /vd de- 
scribes the importance of convection relative to diffusion over 
the vertical distance H. 

First, we consider the finger dynamics from instability on- 
set to steady state (i.e. Fig. [lj> to d). After onset, several 
of the initially shorter fingers stagnate and disappear which 
we interpret as the longer dominant fingers out competing the 
shorter fingers and preferentially draining dense mixture from 
the diffusion interface. Later, the remaining dominant fingers 
move laterally along the interface, collide, and merge result- 
ing in the coarse pattern of Fig. [T}i. When plotted versus the 
scaled depth z/S of the longest finger, S/X shows consistent 
behavior among all the experimental runs (see Fig. [2^). Here, 
8 is the measured finger width. Near z/S ~ 1, i.e. near the 
instability onset, S/X ranges between 0.22 and 0.27 consis- 
tent with the onset value of S/X « 0.27 computed from the 
results in Section 4 of [6|. (Determining the onset A is rela- 
tively simple, however the limitations of our imaging system 
prohibit us from determining the onset time.) For larger z/S, 
8/X decreases as the pattern coarsens and appears to asymp- 
totically approach S/X ~ 0.04. We expect this behavior for 
two reasons. First, if 8/ X — > 0, then the water flux per unit 
interface length would tend to zero for large z or H which is 
unphysical. Second, the coarsening due to finger merging is 
eventually balanced by the nucleation of new fingers as the 
diffusion interface between two now distant fingers thickens 
and becomes unstable (see Fig{T}l), and for large H, this in- 
stability should be governed by local processes independent 
of H. After the longest finger extends to H, we observe some 
fluctuations in 8/X, but the time average is close to the S/X 
near z/S = H/8 in Fig. [2^. We conclude that interactions of 
the fingers with the bottom of the cell do not affect the pattern, 
and 8 is an appropriate length scale for coarsening as opposed 
to H. 

Mixing of the water and PPG is important for understand- 



ing the convective process and assessing its effectiveness to 
securely sequester CO2 in brine reservoirs. We determine the 
amount of water mixed from the experimental images by mea- 
suring the area A between the initial interface and the mid- 
dle of the subsequent diffusive interface (as in Fig. [T]l) as a 
function of t* 03] - Figure^ shows A/LH vs t* for several 
runs with the same K but different H. After an initial slow 
increase (13), A/LH increases linearly with t* indicating a 
steady mass transfer rate. At longer times in runs with small 
H (e.g., the Ra = 8, 700 data), A(t*)/LH saturates due to 
an increase in mixture density in the lower cell which reduces 
Ap and shuts off the convection. 

We express the mass transfer rate in dimensionless form 
via the Nusselt number Nu — m/ p(D/H)bL which repre- 
sents the actual mass flux normalized by the diffusive mass 
flux over cell depth H. Rewriting, Nu = Ra d(A/LH) /dt*, 
i.e. Nu is simply Ra times the slopes from Fig. [2] which are 
plotted versus Ra in Fig. [3] When fit to a power law, Nu 
shows a dependence of Nu = (0.045 ± 0.025)i?a a76±006 
which is similar to results from numerical simulations of ther- 
mal convection in three-dimensional porous media [ 16 1 where 
Nu rj 0M7Ra 9 for 10 3 < Ra < 10 4 . Comparing values 
at Ra — 10, 000, we find Nu w 49 and thermal convection 
simulations give Nu ss 67. Although the comparison with 
simulation is reasonable, we seek to go further by checking 
the internal consistency of the data in Figs. [2^ and [3] and us- 
ing these results to outline a new stability problem that should 
provide a prediction of Nu. 

The mass flux of water through the fingers and into the 
lower cell is given by pc w nSvf where n is the number of fin- 
gers extending from the interface and vt and c w are the ve- 
locity and water concentration of the fluid in the finger cores. 
Scaling this mass flux for comparison with Fig. [3] we find 
Nu/Ra — c w (S/X)(vf/v c ). To estimate Vf/v c , we con- 
sider a set of equally spaced fingers and ignore small density 
differences except in buoyancy terms. The balance of mass 
flux across a horizontal line cutting through the fingers re- 
quires VfS = (A — S)vb, where Vf, is the upward velocity of 
pure PPG between the fingers. Integrating the vertical compo- 
nent of Darcy's law around a closed rectangular contour with 
one vertical leg in a finger and the other in the pure PPG we 
find Vf = v c — (vb/vf)vb- Combining the previous results 
yields Vf/v c = [1 + (vb/vf)5/(X — S)]^ 1 allowing us to ex- 
press Nu/Ra (from above) as a function solely of S/X, i.e. 
Nu/Ra = (S/X)/[l + {v h /v f ){8/X)/{l - 8/X)] = F[8/X] 
(vb/vf is fixed in this experiment). In the quasi steady state 
after the fingers reach z — H, we can instead interpret the 
fitted curve 8/X = G[z/5] in Fig.[2ji as 8/X = G[H/5}. How- 
ever, H/8 = Nu so that 8/X = G[Nu]. Combining, we find 
Nu/Ra = F[G[Nu]] which is numerically inverted to find 
Nu versus Ra as shown in Fig.|3]for c w = 0.3 and 0.5, i.e. 
c w at maximum mixture density and the maximum c w where 
the mixture density is still greater than pppc- 

The rough agreement in Fig. [3] supports the validity of this 
simple model of mass transfer and the relationship between 
Nu and the pattern coarsening (S/X) data in Fig. [2^. How- 
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FIG. 3: (Color online) Nu vs Ra. Different values of H are H = 
1.27 cm (open circles), H = 2.54 cm (filled triangles), and H = 
5.08 cm (filled squares). The solid line is the best power-law fit 
of the form Nu — 0.045i?a 0,76 . The dashed and dotted lines are 
Nu computed from the simple mass flux model described in the text 
using the fit to 8/ A given in Fig 2a. The two curves assume different 
values for the concentration of water in the core of the descending 
fingers. 

ever, a full understanding of the mass transfer must include a 
theoretical model of S/X. Here, we describe a few relevant 
physical process in this model, but leave the detailed deriva- 
tion to future work. In our view, the core of this model in- 
cludes a new stability calculation that predicts the distance 
between two established fingers at which the intervening dif- 
fusion layer becomes unstable to small perturbations. How- 
ever, this layer is already experiencing significant transverse 
flows due to fluid drainage from the layer into the fingers — 
a significant complication not present in the instability of the 
original quiescent diffusion layer. Taylor dispersion due to 
shear along the direction of b in Fig. [T^ will be quite effective 
in smoothing out concentration perturbations. As a result, the 
instability will preferentially nucleate at the symmetry point 
between the fingers (as is seen in Fig. [TJl and in many other 
occurrences of this process) where the flow stagnates and Tay- 
lor dispersion has the least effect. However, a nascent finger 
of finite transverse size at the stagnation point will be spread 
out laterally due to streamwise velocity gradients in the layer. 
Both of these mechanisms behave like an enhanced diffusion 
process which, by analogy with the quiescent diffusion layer, 
will increase the onset wavelength consistent with the small 
8/X observed in quasi steady state in Fig. [2^. 

In conclusion, our measurements of 6/ A at early times are 
consistent with models of this initial instabilityJ3H6l support- 
ing the validity of these models. We have measured both the 
finger coarsening dynamics and steady-state convective mass 
transfer and outlined a new stability problem to predict the 
mass transfer. Our experimental results are useful for un- 
derstanding the impact of reservoir parameters on the effec- 
tiveness of CO2 sequestration. Permeability K accounts for 
much of the variability in Ra between different saline reser- 



voirs. A recent survey of reservoirs in the Alberta basin lfTTll 
found a maximum Ra w 1,000. Extrapolating our best-fit 
power law to this lower Ra yields Nu sa 10. However, other 
reservoirs |fT8ll have significantly higher K with Ra s» 20, 000 
and Nu as 80. Assuming a reservoir depth of 50 m, the char- 
acteristic time for diffusive mixing is td ~ 16, 000 years. At 
Nu s» 80, the convective process observed in this work would 
speed up the mass transfer and reduce this time to approxi- 
mately 200 years. 
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